Exact relaxation in a class of non-equilibrium quantum lattice systems 
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A reasonable physical intuition in the study of interacting quantum systems says that, independent of the 
initial state, the system will tend to equilibrate. In this work we study a setting where relaxation to a steady state 
is exact, namely for the Bose-Hubbard model where the system is quenched from a Mott quantum phase to the 
strong superfluid regime. We find that the evolving state locally relaxes to a steady state with maximum entropy 
constrained by second moments, maximizing the entanglement, to a state which is different from the thermal 
state of the new Hamiltonian. Remarkably, in the infinite system limit this relaxation is true for all large times, 
and no time average is necessary. For large but finite system size we give a time interval for which the system 
locally "looks relaxed" up to a prescribed error. Our argument includes a central limit theorem for harmonic 
systems and exploits the finite speed of sound. Additionally, we show that for all periodic initial configurations, 
reminiscent of charge density waves, the system relaxes locally. We sketch experimentally accessible signatures 
in optical lattices as well as implications for the foundations of quantum statistical mechanics. 

PACS numbers: 



I. INTRODUCTION 

The study of the non-equilibrium properties of quantum 
many body systems has recently entered a renaissance era. 
This has been motivated, in part, by recent experimental de- 
velopments; the rapid progress of experiments involving ultra- 
cold atoms in optical lattices, with their high degree of con- 
trol and long coherence times, has opened the door to precise 
experimental studies of the dynamics of strongly interacting 
quantum systems |Q]]. In addition, questions of relaxation and 
thermalisation for non-equilibrium systems are again receiv- 
ing attention from foundational perspectives. This is partly 
triggered by intuition from quantum information theory where 
maximally or almost maximally entangled states emerge from 
appropriate distributions of random states [2-4]. 

One particularly fascinating setting which has recently re- 
ceived intensive study is that of quenching, that is, a sudden 
change of interaction strength [5-12]. Several explanations 
for, and numerical studies of, quenched dynamics have grad- 
ually led to the formulation of a rather general body of theory 
and conjectures; it has been mooted that if the system starts 
in the ground state of one Hamiltonian then certain properties 
such as correlators of the system should relax to an analogue 
of the thermal state of the new Hamiltonian after a quench 
[5, 7-9]. Thus we are motivated by these observations to for- 
mulate a conjecture. This local relaxation conjecture states 
that a system should locally relax to a steady state, respecting 
conserved quantities of motion. 

The local relaxation conjecture may sound suspiciously like 
a violation of unitarity as the global state must, of course, re- 
main pure throughout the dynamics. However, a reasonable 
physical intuition which explains why there is no contradic- 
tion is that for a small block of sites the rest of a system acts 
like a reservoir and thus allows the site to maximise its entropy 
subject to the constraint that the energy is preserved. Indeed, 
the full explanation for the emergence of a local steady state 
during the course of the quench may be intuitively described 
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FIG. 1: Intuitive picture of the relaxation process in the quenched 
Bose-Hubbard model as discussed in the main text: For any lattice 
site i (or any block of sites) within a cone (dark grey) defined by the 
speed of sound excitations significantly contribute to the local mix- 
ing of the state at the site. Contributions from outside this cone (light 
grey) are exponentially suppressed. The incommensurate influence 
of the lattice sites in the cone gives rise to a relaxation to maximal en- 
tropy for large times t. When the excitations have travelled through 
the entire lattice, recurrences will occur. 



along these lines: as time evolves the system becomes corre- 
lated ifToll - from each site a wavefront moving at the speed 
of sound for the lattice emerges carrying information. As 
time progresses more and more excitations will have passed 
through a given site, see Fig. 1 . The cumulative effect of these 
successive excitations results in an effective averaging pro- 
cess; information stored in one site becomes infinitely diluted 
across the lattice as time progresses. A similar intuition has 
been recently emphasised by Calabrese and Cardy lfl3ll in the 
context of quenching to a critical system. 

In this work we introduce a physical setting where the local 
relaxation conjecture can be studied analytically. We find an 
exact convergence of all local states to a steady state for long 
times in the quenched Bose-Hubbard model. We imagine the 
system starts in a Mott insulator phase and is then suddenly 
switched to the deep superfluid regime. Our method is self- 
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contained and physically motivated, and is valid also for finite 
system sizes relevant to experimental settings. In the course 
of the proof we are also able to quantify the dilution of infor- 
mation throughout the lattice. Thus, our approach is consid- 
erably simpler and more physically intuitive than a potential 
approach based on the C* -algebraic arguments developed in 
Ref. Hill to study the local relaxation of free fermions and 
bosons freely moving in R™ and classical systems lfl5ll . Inter- 
estingly, the convergence is not only true in the time average, 
but actually true for any large instant of time. This is in con- 
trast to the recent approaches developed in quantum informa- 
tion where, in order to study this problem, it would seem 
necessary to consider the time-averaged local state rather than 
the local state itself. Additionally, our approach does not re- 
quire the system being quenched to to be critical — we ob- 
tain the same results regardless of whether the system being 
quenched to is critical or not. For blocks of sites we also find a 
relaxation, but not to the thermal state of the new Hamiltonian. 

II. QUENCHED BOSE-HUBBARD DYNAMICS 

We start from a Bose-Hubbard Hamiltonian, modeling, for 
example, a dilute gas of ultra-cold atoms which are Bose- 
condensed in an optical lattice lfl6ll . In one dimension (gen- 
eralizations of our findings to higher dimensions are entirely 
straightforward) it reads H = Hs F + U, where 

N N 

H SP = -jJ2blbj, U= -^n^-l)-^^- W 

i=l -t=l 

Here b\, .. . ,bjsi denote the bosonic annihilation operators, 
and hi = SjSj. The coefficient J is the hopping matrix el- 
ement, governing the strength of hopping to neighboring sites 
((i, j) indicates summation over nearest neighbors) and U de- 
fines the strength of the on-site interaction. Finally, the chem- 
ical potential p controls the particle number. For theoretical 
simplicity we shall assume that that the system is translation- 
ally invariant with periodic boundary conditions, so that the 
underlying lattice is a ring. This Hamiltonian exhibits two 
distinct phases as J and U are varied. When the hopping dom- 
inates, that is J ^> U, the system is in a superfluid phase. For 
a dominant on-site repulsion the ground state is a Mott insu- 
lator. 

We imagine that the system is initially held at chemical 
pote ntial p and is in the ground state in the Mott regime 
111611 . In this deep Mott phase (corresponding to J = 0) 
the ground state vector is a product \m)® N of individual 
number state vectors \m) of to bosons at each site. Here, 
m = max{0, [{n/U + l)/2]}, where [•] denotes the closest 
integer to the value in brackets. We then imagine that the sys- 
tem is rapidly quenched, at t = 0, to the strong superfluid 
regime J ^> U. We model this quench by an instantaneous 
change in the Hamiltonian to the regime U = 0. Thus the 
system state vector at time t, \ijj(N,t)), is given by (we set 

h= l) 

\1i(N,t)) = e- u6 "\m) 9N . (2) 



III. RELAXATION FOR A SINGLE SITE 

As outlined earlier, we are interested in whether subsystems 
consisting of single sites or blocks thereof equilibrate, and, if 
so, in what sense. For clarity we start by proving relaxation 
for the quantum state Qi (N, t) of a single site i in the lattice 
L = [1, . . . , N], and then extend this result to multiple sites 
and different initial conditions. The state of site i is given by 
a partial trace 

Qi(N,t) = tr L \ {i y[\ip(N,t)){ip(N,t)\]. (3) 

Throughout this paper we will describe the system in phase 
space, making use of the characteristic function to represent 
the state Qi. It is defined as 

Xi(a;N,t) = ix\§ i (N,t)e c& *- a ' tt ], (4) 

where a £ C. We now turn to showing that we indeed lo- 
cally, at each site i, find convergence to a state that maximises 
the entropy: We prove that the state Qi(N, t) approximates a 
Gaussian state qg arbitrarily well. We will find that the first 
moments of qg vanish, tr[^c^i] = 0, and its second moments 
are all conserved and identical to the initial ones, 

tv[g G bibi] = (to | | to) = 0, (5) 
tT[g G hi] = {m\fn\m) = m. (6) 

Gaussian states maximise the local entropy for fixed second 
moments, which are constants of motion for the initial states 
\m)® N under consideration. For more general initial states 
exhibiting time dependent second moments see Section lVTl 
Our main result can thus be summarised as 

ft(JV,t)-fg G , (7) 

where the limit of large N is taken first and then the limit 
of large t is taken — note that there is no time average in- 
volved. More precisely, for any e > and any recurrence 
time ^Recurrence > we can find a system size N such that 

\\Bi(N,t) - Q G \\v < e forte [t 

Relax: ^Relax ^ Recurrence ]• (8) 

Here, the relaxation time tR e i a x is governed by the hopping 
strength J, which defines the speed of sound of the system, 
see Eq. d24l ) and Fig. Q] For finite N recurrences occur for 
times larger than <R e cun-ence, however these can be shifted to in- 
finity for large N. Closeness of Qi(N, t) and qg is measured 
in trace norm which of course means that all local expectation 
values are the same as for the relaxed state. Another way of 
interpreting the relaxation to a Gaussian state is that the entan- 
glement between the site i and the rest of the chain becomes 
maximal. 

For clarity, we deliberately discuss this case of a single lat- 
tice site first in great detail. We will show Eq. ([8]) by evaluat- 
ing the corresponding limits for the characteristic function \i, 
and showing that it tends to a Gaussian in a for large times t 
and lattice sizes N. Pointwise convergence in the characteris- 
tic function has been shown to imply trace norm convergence 
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of the corresponding density operators 1124 12511 . so this is suf- 
ficient to prove Eq. {8). There are four main steps in the ap- 
proach: (i) We take and expand the logarithm of the character- 
isic function, then (ii) we identify certain terms quadratic in a. 
(iii) The magnitude of the remaining terms are then bounded, 
and finally (iv) we make use of these bounds to evaluate the 
limiting behavior of the characteristic function. Readers inter- 
ested in the final result may safely skip this part of the argu- 
ment, which is completed by Eq. ( f32b . 

Employing the cyclic rule of trace, we find 

Xi(a;N,t) = tr[ C - i ^ SF |m)® JV (m| 8Ar e i ^ SF e^ t - Q * S - 

= { m \ 9N fP A '*e ai t- a " 6 *e- it6 *\m) 9N . (9) 

Now, using the Baker-Hausdorff identity (or, alternatively, by 
solving Heisenberg's equation of motion) one finds for the 
Heisenberg representation of bi 

N 

e itH SF g. e -im SF = J- [V(N,t)] itj bj, V(N,t) = e~ itJ ( N \ 

3=1 

(10) 

where the N x N circulant matrix J{N) represents the hop- 



ping operator Hsf- Its entries read 



-J8, 



dist(i ,j 



i, where 



we denote by dist(i, j) the distance between two sites i and j. 
Due to this circulant structure J{N) may be diagonalized via 
a discrete Fourier transform to give eigenvalues 



\ k = -2Jcos{2irk/N) , 



and thus 



1 N 



(ID 



(12) 



k=l 



The characteristic function in Eq. (0 is then given by 
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I at + \ TT / I a\V(N,t)]* M-a'lVlNMi ,t>,- 1 \ 

Xi( a ! N ,t) = II \tn\c 1 v ' " ,J 3 \m). 

3=1 



13) 



A straightforward calculation (see, for example, Ref. IU8I1 ) 
shows that the expectations appearing on the right-hand side 
evaluate to Laguerre polynomials and we may write 
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Xl (a;N 7 t) = e-l a l 2 / 2 JjL m (|A,,-(iV,t)| 2 ), (14) 

3=1 

with [3ij(N,t) = a[V(N,t)]*j (we will often assume im- 
plicit dependence on t and N and simply write (3ij and Vij). 

The matrix elements Vij encode all the information about 
the dynamics. We identify each row of this matrix as the quan- 
tum state of a free particle initially localised at site i which is 
hopping on the ring. As time progresses the initially localised 
particle disperses rapidly throughout the ring. In deriving the 
bounds of step (iii) it will be essential that we can identify for 
any e > those times t and lattice sizes N for which we can 
bound \Vij\ < e. This part has, in turn, two main ingredi- 
ents: We make use of (iiia) a "central limit type argument" 



as well as of (iiib) a "Lieb-Robinson type argument". In this 
way, we also have a very clear handle on the case of finite N, 
which is the situation encountered in any experiment. These 
two regimes of (iiia) and (iiib) correspond to lattice sites close 
and far away from the distinguished site i: For sites labeled 
j with dist(t, j) < L(t) for some L{t) we are able to grasp 
the dynamics using a central limit argument. These are the 
sites actually contributing to the relaxation process (see Fig. 
1). Sites with dist(i, j) > L(t) also contribute to the local 
dynamics at site i, but this contribution is exponentially sup- 
pressed due to the finite speed of sound for the system. 

We start by investigating the sites with negligible contribu- 
tion, case (iiib). As Jlj = for dist(z , j) > 1, we have that 
[J n ]i.j = for dist(z, j) > n and thus 



E 

y. Hi) 



(-it)* 



Vij = \e-* Jt ]. . 

n=0 



(15) 



[J n ] 



n>dist(i,j) 



Now, for any matrix M one has \Mij\ < \\M\\, where 
indicates the operator norm, and ||,7|| = 2J. Hence, 



n! ' — ; \ n 

n>dist(ij') n>dist(i,j) 



< 



E 



n>dist(i,j) 



6Jt 



dist(i,j')y 



(16) 



where we used the fact that n\ > (n/3) n . We thus finally find 
that contributions from sites j with dist(i, j) = d > 6Jt are 
exponentially suppressed: 



h3 



< 



{6Jt/d) d 
1 - QJt/d 1 



(17) 



independent of the system size N . This is the intuition pro- 
vided by what are known as Lieb-Robinson bounds: Sites be- 
yond the cone defined by the speed of sound will not signif- 
icantly alter the state at site i (see again Fig. Q]). From this 
expression we see that, if we require \ Vij \ < e, then it is suf- 
ficient that dist(i, j) = d > L(t), where L(t) is given by the 
solution to (6Jt/L(t)) L ^ = e(l-6Jt/L(t)). A crude bound 
on L(t) may be obtained by noting that 



{6Jt/d) d < QJt/d 



1 - QJt/d ~ 1 - 6Jt/d' 
i.e., for given e > Owe have 

\Vij\ < eforalH, j : dist(i,j) > L(t) = 6Jt 



1 + e 



(18) 



(19) 



We now turn to the case where dist(i, j) < L(t), case (iiia). 
Writing Vi.j = V\ for / = i — j, we recall that 



1 N 



3 2it J cos(27rfe /N) 2-rrikl / N 



(20) 



fe=i 
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FIG. 2: The entries of the matrix governing the dynamics of the 
system are arbitrarily small — the system is locally arbitrarily close 
to a Gaussian state — for sufficiently large TV and a time interval 
[tReiax, iRciax + ^Recurrence] (compare Eq. <(24j> and Fig. [TJ. After ex- 
citations have traveled through the whole lattice recurrences occur at 
t = tRecurrence, which is given by the speed of sound. For TV — > do no 
recurrences occur, tRecurrence — * <x>, and locally the system relaxes to 
a Gaussian for large times. 



In the limit TV — > oo this approaches an integral representation 
of the Bessel function 



e i»cos(</)) e i/0 



(21) 



up to a phase r . For finite N, the Vi can thus be thought of as 
Riemann sum approximations to i' J;. Our strategy for choos- 
ing finite TV, t such that \Vi\ < e will make use of a bound 
on the error involved in such an approximation, together with 
the bound \Ji(x)\ < x' 1 ^ for all x > and all I JH]. Con- 
sider the quantity \Vi — i Ji(2Ji)\. On using the Riemann sum 
approximation, we find 



IVJ - i l Ji{2Jt)\ < 2ix 2 {2Jt + l)/N, 



(22) 



and by combining this with the bound on \ Ji(2Jt)\ we obtain 

\Vi\ < ^-(2Jt + L(t)) + (2Jt)- 1 / 3 , (23) 

complementing the bound in Eq. (T% . 

Combining both bounds we find that for all i,j (see Fig. [2] 
where we plot maxij | V^^- 1 versus TV and t) 



4 TV e 2 

\V itj \ < e for aUt : - < Jt < — — . (24) 
e d 47r z 8e + 6 



For any e > 0, this provides sufficient conditions for satisfy- 
ing \Vij\ < e for all i, j and finite N and t. This ends step 
(iii), the intuition of which is summarized in Fig.Q] 

We now turn to our original aim, that of showing that the 
characteristic function Xi tends to a Gaussian in a. The steps 
in the proof follow the basic structure of the classical proof 
of the central limit theorem lfl7ll . We start by taking the loga- 
rithm of the right hand side of Eq. (TBl i. Expanding the loga- 
rithm, we then make use of the observations on \ Vij | obtained 



above, together with certain properties of Laguerre polynomi- 
als. 

We can always find appropriate t and N such that < 
L m {\Pi,j\ 2 ) < 1 as j | can be made arbitrarily small for 
given a and all see Eq. d24b . Thus we may write 

log n mi Aj i 2 ) = - E E ( r )j ■ (25) 

3 = 1 3=1 fe=l 

The Laguerre polynomials are defined as 



L m{x) = E 



71 = 



m 
m — n 



l — mx + l m (x), (26) 



where we defined l m (x) for later use. Hence 

N N N 

l0 g n^m(IA,3-| 2 ) = -^El^f + E^d^l 2 ) 
3=1 3=1 3=1 

_ ^^ (i~L„m,\ 2 ))\ (27) 

3=1 k=2 

From the definition of and the unitarity of V, the first 
summation evaluates to — m|a| 2 , which is quadratic as de- 
sired, and the remainder of the proof consists of showing that 
the remaining terms tend to zero for given a. 

Since \ = \aVi.j |, we may use the results above to en- 
sure that for any fixed a the \j3i j \ are as small as desired, and 
in particular we shall first suppose that | < 1. In this case 
we have |^m(|Aj| 2 )l < | 4 ^m(-l), and it follows that 



JV 

^ ] ^m(|A,3 
3 = 1 



N 

<l m (-l)\a\±J2\Vi,j\ 4 

3=1 

= :l m (-l)\a\ A V. 



(28) 



The large TV behavior of this term will evidently be deter- 
mined by the corresponding behavior of V, for which we al- 
ready know that it tends to zero if we first let TV and then t go 
to infinity, see Eq. (l24l i. 

Now we study the last term in Eq. ( f27T >: 



E 

fc=2 



(i-L m m,j\ 2 )) k 



(29) 



fc=i 



k + 1 



< 



- {1 - L m (\(3, u \ 2 )) logL m (\p itj \ 2 ). 



By again expanding the Laguerre polynomial as above we 
have 

-logL m (|/3, M | 2 ) < -log(l-m|A,3f-|«m(|A,3-| 2 )l) 
< -log(l- |a| 2 ^(m + / m (-l))) 
=: A, (30) 
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which is a sensible bound if we require that |a| 2 y(m + 
l m (— 1)) < 1, which can again be ensured with appropriate 
N and t. Furthermore, V — > 0, i.e., also A — > 0. We thus 
finally obtain the upper bound 

N 

AJ2 (i - Ml/M 3 )) = AEf =1 (HA,il 2 - imdA.il 2 )) 

i=l 

< A(m|a| 2 + |a| 4 VZ m (-l)) (31) 

for the absolute value of the last term in Eq. (1271) . 

Equipped with these bounds, we arrive at the desired state- 
ment: We find for any aeC, 

N, t) = e -(™+V2)l«l 2 + ^ t)j (32) 

where \f(a;N,t)\ can be made arbitrarily small for all t € 
[£Reiax,*Recun-ence + *Reiax] and a suitable choice of the lattice 
size N, see Eq. d24l i. Since pointwise closeness for the char- 
acteristic functions of two quantum states translates to close- 
ness in trace norm for the states themselves lf24ll25ll . we have 
completed the argument for a single site. 

We therefore find that the state of a single site, in the limit of 
large t, maximises the entropy for the given initial second mo- 
ments: Gaussian states have this extremal property of having 
maximal entropy. Another way of saying this is that, subject 
to this constraint, it is maximally entangled with the rest of the 
chain. Since we can bound the time scale on which relaxation 
occurs - the single site relaxes to a Gaussian state at the rate 
(2Jt)~ 1 ' 3 - our argument is robust under small perturbations. 

This result also gives a clear handle on the situation for fi- 
nite system sizes: For any finite system size N, there will of 
course eventually be recurrences. The time scale for such re- 
currences to occur in the (JV, t) plane is governed by the speed 
of sound and the system size N: when t < cN, where c is the 
speed of sound, there will be no recurrences, and our argu- 
ment will apply, giving rise to the time ^Recurrence, see Eq. d24i >. 
Thus we also arrive at a precise way to describe relaxation 
for infinite system size: locally the system "looks exactly re- 
laxed" for large times, apparently giving the impression of an 
"irreversible dynamics". The same result also holds for cor- 
relation functions. This local relaxation is perfectly consis- 
tent with the unitary dynamics of the full system; the whole 
system never truly relaxes to an equilibrium state as it must 
always preserve the full memory of the initial condition. 

The Gaussian character of the state can also be viewed as 
resulting from a different type of central limit theorem (cf. 
Ref. 02311 ). The intuition is that time evolution can be thought 
of as a sequence of of beam splitters acting in parallel, i.e. a 
quantum quincunx or Galton's board (this follows from the 
fact that the dynamics are well-approximated by a quantum 
cellular automaton). After the bosons bounce through the 
quincunx we should arrive at a Gaussian characteristic func- 
tion in phase space describing the state of any single site. 



IV. RELAXATION OF BLOCKS 

The equilibration for multiple sites, which includes blocks 
of lattice sites, may be shown using identical techniques, and 



here we briefly describe the necessary, but minor, modifi- 
cations. The characteristic function as a function of a = 
(ai, . . . , a s ) T £ <C S for a set of sites S = [1, . . . , s] with 
reduced state q$ is defined as 



Xs(a) = tr Q S {N,t)e 



IP E? =I M-c 



(33) 



After the quench the time-evolving bosonic operators are 
given by Eq. (1101 . and the characteristic function becomes 



- p-aWa 



TV 



Y[L m (\PsAN,t)\ 2 ), (34) 

3=1 



where the 0s j are now given by aiVy(JV,t), Note 

that otherwise Eq. (fl4l is unchanged, and the argument can 
be followed as before to show pointwise convergence in a. 
Using arguments essentially identical to the ones in the proof 
for a single site, we find in the same sense a convergence 



Q S (N,t)->Q%' 



(35) 



for large times t and large system sizes N. That is, the block 
tends to a product of maximal entropy states. 



V. COMPARISON WITH THERMAL STATES 

Let us contrast the state arising from our dy- 
namical relaxation with the thermal Gibbs state 

e -(j^JV)/T/ tr [ e -(ff SF -^iV)/T] (setting kg = 1} of 

the new Hamiltonian Hgp: for a thermal state at temperature 
T > 0, we find the correlations 12211 



(e 



(J-Hl)/T 



iy 



(36) 



The thermal state of Hsf is always correlated for all temper- 
atures T > and all chemical potentials fx fixing the particle 
number. In other words, although the system does mix and we 
do arrive at a state with local maximal entropy, it does not re- 
lax to the thermal state of the new Hamiltonian: the conserved 
second moments prevent this from happening. The relaxation 
should be thought of as resulting from a mixing process rather 
than of a physical thermalisation process. Of course the sec- 
ond moments are locally preserved by the dynamics, so under 
this constraint the entropy is indeed maximised. 



VI. MORE GENERAL INITIAL STATES 

In fact, an argument along the lines of the above proof 
can be carried out whenever the the initial state has suffi- 
ciently rapidly decaying correlations, reminiscent of the sit- 
uation in the classical ca se 1 1511 . Details of this argument will 
be presented elsewhere 12611 . For the purposes of this article, 
we will now consider more general initial states of the form 



®iL\\ m i)' allowing for different particle numbers at differ- 
ent sites - allowing, e.g., for configurations occurring in har- 
monic traps and periodic distributions with periodicity differ- 
ent from one (realized, e.g., by superimposing a second lat- 
tice with different wavelength), reminding of charge density 
wave or checkerboard phases. For these inhomogeneous ini- 
tial states the second moments are no longer preserved, but for 
certain initial states we can nevertheless prove a convergence 
of second moments. While convergence of second moments 
is necessary for a local relaxation, it is certainly not sufficient. 
However, as before, we show that local states relax to a Gaus- 
sian state - we show that the characteristic function of the state 
tends to a Gaussian in a as higher cumulants can be shown to 
vanish just like in the homogeneous case. For initial states 
0)) = <%>iL-i\mi), the expression for the characteristic 
function is 

JV 

Xs(t*;N,t) =c- a ' a / 2 l[L mj (\(3s. J (N,t)\ 2 ), (37) 

3=1 



with \(3s.j \ as above. Eqs ( 128H29t and the surrounding discus- 
sion must then be modified to incorporate the different rrij, 
for example in the quantities l m (— 1). All such quantities can 
be appropriately modified by taking their supremum over all 
the different rrij. We then find pointwise convergence of the 
characteristic function to 



t— >oo JV— >oo 



lim lim xs(a;A r ,i)=c jQ/2 

T -oc 

N 

exp(— lim lim / WjlAs 1 1 2 ) , (38) 

V t— >oo N— too ' / 



3=1 



showing that there can be initial configurations for which the 



state never relaxes. For initial states with m>. 



m for almost 



all j, we always find relaxation to a steady state as we also do 
for P-periodic initial configurations: Consider an initial state 
with rrii = rrii-p, completely determined by mi, . . . , mp. 
Then, for a block S of s sites, 




-200 
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FIG. 3: Second moments (b\b 1+ i) (left) and {Sg^+i) (right), 
I — 0,1,5,10,15,20, for an initial periodic configuration 
|0, 1, 0, 1, • ■ ■ ,0, 1) on a lattice of size N = 300. For larger 
times recurrences occur (not shown), while for N — > oo, one has 

lim t ^ 00 (k-&j) = Sij/2. 



which yields 



N/P-l 

E 

3=0 



Et/* t/ ^ k > 1 i 1 (z+fc)+^ 

v k ,p j+p v lt p j+p = — + - ^ c [p[ 4 



p-i 



z=-P+l 
z#0 



x J_ -2itsin(-7rz/P)V 2nir(l-k)/N 



(42 



Identifying the last line with [V(N,—2t8in(irz/P))]i t h, 
which tends to zero as t and N become large (see the proof 
for the relaxation of a single site), we find also in this scenario 
that the characteristic function tends to a Gaussian: 

lim lim X s{a;N,t) = B -W3+*>)<*\ h ^ 

t^oo N—too 

(43) 

where m = ^2 p=1 m p /P. While initially different sites are 
uncorrelated, correlations build up over time and finally go to 
zero for large times (see Fig. [3j and the local state relaxes to a 
direct product of the same Gaussian at each site. 



JV s P N/P-l 

^Z m 0s,j\ 2 = E E m P E V k,Pj+pVi,Pj+ P - 

3=1 M=l P=l 3=0 

(39) 

Now, 



N/P-l 

E 

3=0 



/P-l 1 N 

E v * \r _ „it(A r -A 3 )„2iri(W 

Vk,Pj+pVl,Pj+P ~ T72 e e 



jV 2 

r.s — 1 
N/P-l 



E e 

3=0 



2mPj(r-s)/N 



sk)/N 



(40) 



and the last terms evaluate to a comb of Kronecker deltas, 



N/P-l 

E 

3=0 



^iPj(r-s)/N 



N 

p 



p-i 



E *r- 



s,zJV/Pj 



(41) 



z=-P+l 



VII. INTEGRABLE SPIN MODELS AND FERMIONS 

To further illustrate the generality of our approach, we 
sketch the situation in other integrable models. We consider a 
XY spin chain, 

H = - 1 -y( -±-^ + i— - - y at, 

(44) 

again in the setting of a quenched interaction. By virtue of 
the familiar Jordan- Wigner transformation, the Hamiltonian 
is equivalent to a system of spinless fermions: 

" ^E [fl^jh - umJ] + /,//,,./• - /;//,,/;) , 

i,3 

(45) 
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where An = A, Ai 



-1/2 for dist(i,j) = 1, B t 



—Bj i = 7/2 for dist(i, j) = 1 and zero otherwise. 

We start with the system initially in one of the states 
\m)® N , m = 0,1, corresponding to the magnetized spin 
states 1 1)® N , I T)®"^, an d study time evolution under the 
above Hamiltonian. As the initial states are fermionic Gaus- 
sian states, the system stays Gaussian for all times t. Hence, 
we only need to consider seconds moments, which are most 
conveniently written in terms of Majorana operators, m, = 
(fi + fi )/v / 2, m N+l = i(fl - fi)/V2. Then one may collect 
second moments in a correlation matrix 7 (AT, t) with entries 
[j(N,t)]ij — u[g(N,t)rhithj]. Time evolution of the initial 
state \m)® N under the above Hamiltonian then amounts to 



-y(N,t) = e tjH 7(A,0)c 



-tH 



H 



V 
-V T 



where V = A + B with eigenvalues 



Afc(V) = A - cos(27rfc/A) - i7sin(27rfc/A). 



(46) 



(47) 



A tedious but elementary calculation delivers the resulting 
second moments 



where 



1 + i(-l) r 



Mi(JV,t) M 2 (N,t) 
-M%{N,t) Ml 



M x (N,t) = sin(2|V|i) 



B 

\V\' 



M 2 (N,t) 



+ sm 2 {\V\t) 



W\ 



1 



(48) 

(49) 
(50) 



from which we immediately see that for the isotropic XY 
model (7 = 0) second moments are a constant of motion, 
7 (AT, t) = 7(AT, 0). We assume 7 / from now on. Ob- 
viously, as in the bosonic case, the total state stays pure for 
all times and its entropy is thus always zero, however, de- 
pending on the system parameters (it will turn out that it 
depends on whether the system being quenched to is crit- 
ical or not - in contrast to the bosonic case, where relax- 
ation occurs independent of the system being critical or not), 
we find exact relaxation of second moments. To this end, 
consider the second moments of the state of a single site, 
Qi(N,t) — trjA/ji [@(N, t)]. It is completely determined by 
the reduction of the correlation matrix to site i, i.e., we need 
to consider [Mi(N,t)]ij and [M 2 (N, As the matrices 

A and B commute and M\(N,t) is the product of a sym- 
metric and an antisymmetric matrix, we immediately find 
[Mi(AT,i)] M = = [Mi(A,0)] M . For the sake of clarity, 
we now focus on the Ising model, 7 = 1. It should be clear 
how the following arguments may be extended to explore the 
full parameter space. We find 



[M 2 (N,t) - 1]. 4 = i^sin 2 (|A fe |i) - 1 



fc=i 

N 



l£sin 2 (|A fe |^ 



|A fc | 2 

e 4irifc/JV _ ^ 



(51) 



fc=i 



X c 27Tik/N ■ 



Let us illustrate this for one particular non-critical case, 
namely A = 0, for which one has Ak(V) = —e 27Tlk / N and 
hence 



[M 2 (A,i)]^ = l-sin 2 (t), 



(52) 



which never relaxes. 

It turns out that for the non-critical case, |A| 7^ 1, the sys- 
tem never locally relaxes. This is, of course, no surprise: in 
the case where A = the Hamiltonian we quench to contains 
only commuting terms. Physically, one of three different pro- 
cesses can take place. In the first case, A < 1, the system is 
quenched to a state in the same phase. In this case we do not 
expect local relaxation as the density of quasiparticles of, eg., 
||)® m with respect to the new Hamiltonian is too low: there 
are not enough quasiparticles to lead to an averaging effect. 
In the critical case A = 1, the quasiparticles are delocalised 
and we expect that the density of quasiparticles in | j)® m is 
high enough to lead to a local relaxation. Finally, in the case 
A > 1, the quasiparticles of the new system are strongly lo- 
calised, and appear as pairs in the state | j)® m which then 
oscillate locally. 

Let us now turn to the critical case A = 1, where we expect 
relaxation. We find 

1 N 

[M 2 (N,t) - 1], M = --^sin 2 (2tsin(7rfc/A0), 



k=l 



1 N 

- — ^sin 2 (2tsin(7rfc/A))( 



2mk/N 



(53) 



k=l 



where the first term approaches —1/2 and the second term 
goes to zero for large N, t, i.e., indeed, at each site the system 
relaxes when evolved in time under this critical Hamiltonian. 



VIII. SUMMARY AND DISCUSSION 

In this work, we have investigated the relaxation of Bose- 
Hubbard type systems and related models following a sudden 
quench. The intuition we developed was that within a causal 
cone of neighboring lattice sites, the incommensurate influ- 
ence of the propagating quasiparticles will give rise to a re- 
laxation dynamics. We have considered a setting in which we 
could prove an exact local relaxation to a maximal entropy 
state, providing a guideline of what is expected to happen in 
similar cases. Also, a precise understanding follows for the 
scaling of this effect in the system size, as well as of relaxation 
times. Indeed, the system "looks" perfectly relaxed, while at 
the same time keeping a perfect memory of the initial situa- 
tion. We hence can clarify the difference between an apparent 
- which does happen - and an actual global - which does not 
happen (!) - relaxation. The argument presented here is ex- 
pected to be valid whenever the initial condition is sufficiently 
clustering, so if one has an appropriate decay of correlations. 

From the perspective of kinematical approaches to quantum 
statistical problems, one can interpret our results as follows: 
When randomizing over all possible pure states flU, the local 
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state will have maximal entropy for large systems. Here, we 
arrive at the same result, but not by a kinematical argument, 
but via a dynamical argument, where the mixing is achieved 
through local physical dynamics in a lattice system. It would 
be interesting to further try to combine the kinematical and 
dynamical pictures for more general Hamiltonians. 

Settings similar to those discussed here are readily acces- 
sible in experiments using ultracold bosonic atoms in optical 
lattices Hf], and this is what motivates this work. The general 
guideline is that local or correlation functions will relax. Most 
accessible would be a situation where only second moments 
would have to be monitored, starting, e.g., with a checker- 
board type situation of alternating occupation numbers in a 
Mott state (cp. Section [VT] and Fig. [3}, prepared using two 
distinct optical lattices in a spatial dimension. The idealized 
case of a vanishing interaction U gives rise to a guideline for 
what to expect in realistic settings and also to a general prin- 



ciple: the local relaxation conjecture. It is our hope that the 
present work fosters further experiments along these lines. 
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